// -*- coding: utf-8 -*-

#include <iostream>
#include <fstream>
#include <boost/bind.hpp>

#include "hybrid_model.h"
#include "hybrid_solver.h"

USING_ELISA

VnodePtr HybridSolver::m_odeSolver;

#ifdef HYBRID_PROFILE
long elisa::numSplitTime(0);
long elisa::numHybridReduction(0);
long elisa::numNarrow(0);
long elisa::numNewton(0);
#endif

//--------------------------------------------------------------------------
// class HybridMakeNeighborhood
//--------------------------------------------------------------------------

/**
 */
bool HybridMakeNeighborhood::operator()()
{
	VarImplPtr vip = m_sc->get();

	if (vip->getType() != "TimeImpl") {
		vip->split();

	} else { // TimeImpl
		// (約) 1/2のステップ幅でODE求解を行う。
		// ステップ幅でTimeImplを分割する (3つ以上のステップになった場合、
		// 最初の1つとそれ以後の和をとったものとに分割する)。
		// 左右のIV, ABもTimeImpl内に保存する。

		TimeImpl & ti(*(dynamic_cast<TimeImpl *>(vip.get())));
		Time & time(*(dynamic_cast<Time *>(vip->getSymbolicRep())));

		// TODO: excessの扱いが適当
		VNInterval initT(time.getInitTime() + time.getExcess());
		//VNInterval initT = VNInterval(ti.domain().inf() + time.getExcess());
		VNInterval endT0(ti.domain().inf() + v_bias::sup(time.getExcess()));
		VNInterval endT (ti.domain().sup());

		VNIntVec iv0(time.getInitValue());

		VnodePtr m_odeSolver(HybridSolver::getOdeSolver());

		// Tの下限までIVP-ODEを解く。
		if (v_bias::sup(initT) < v_bias::inf(endT0)) {

			// VNODEの設定。
			m_odeSolver->setHmax(v_bias::sup(endT0) - v_bias::inf(initT));
			m_odeSolver->setFirstEntry();
			m_odeSolver->setOneStep(vnodelp::off);

			m_odeSolver->integrate(initT, iv0, endT0);
		}

		// VNODEの設定。
		// TODO: endTに足りなくないように若干増やしておく。
		m_odeSolver->setHmax((v_bias::sup(endT) - v_bias::inf(initT)) / 2 * 1.00001);
		//m_odeSolver->setHmin(time.getExcess().width() * 2);
		//m_odeSolver->setHmin(1.0e-10);
		m_odeSolver->setFirstEntry();
		m_odeSolver->setOneStep(vnodelp::on);

		//const VNIntVec & iv0 = ti.domain().getIV();
		//const VNIntVec & ae0 = ti.domain().getAE();
		VNIntVec iv(iv0);
		VNIntVec oldIV;

		Interval pivot;
		VNIntVec  ivL, ivR;
		//VNIntVec  aeL, aeR;

		m_cc.initIntegrate();
		bool first = true;

#ifdef _DEBUG
static std::ofstream fot("tight.dat");
// TODO
//printAsRect(fot, initT, iv[ContVar::logIndex]);
//printAsRect(fot, iv[0], iv[2]);
//printAsRect(fot, iv[0], iv[2]);
#endif

		/*
		// 0621
		// count buffered steps
		int sbc(0);
		*/

		while (v_bias::sup(initT) < v_bias::inf(endT)) { // initTとendTはdisjoint

			oldIV.assign(iv.begin(), iv.end());

			m_odeSolver->integrate(initT, iv, endT);

#ifdef _DEBUG
// TODO
//printAsRect(fot, initT, iv[ContVar::logIndex]);
//printAsRect(fot, iv[0], iv[2]);
//printAsRect(fot, iv[0], iv[2]);
#endif

			VNIntVec  ae = m_odeSolver->getAprioriEncl();
			v_bias::interval ti = m_odeSolver->getT();
			
#ifdef _DEBUG
static std::ofstream foi("interval.dat");
// TODO
//printAsRect(foi, ti, ae[ContVar::logIndex]);
//printAsRect(foi, ae[0], ae[2]);
//printAsRect(foi, ae[0], ae[2]);
#endif

			if (first) {
				first = false;
				ivL = oldIV;
				//aeL = ae;
				pivot = Interval(initT);

				ivR = iv;

				break;

			} else {
				/*
				if (ivR.size() == 0) { // 2つ目
					ivR = oldIV;
					aeR = ae;
				
				} else { // 3つ目以降

					// 0621
					if (ae0[0].right() != GAOL_INFINITY && sbc >= SBMax) {
						aeR = ae0;
						break;
					}
					++sbc;

					transform(aeR.begin(), aeR.end(), ae.begin(), aeR.begin(),
						boost::bind(VNIntervalTraits::hull, _1, _2));
				}
				*/
			}
		}

#ifdef HYBRID_PROFILE
		++numSplitTime;
#endif

		//ti.split(pivot, ivL, ivR, aeL, aeR);
		ti.split(pivot, ivL, ivR);
	}

	return true;
}


//--------------------------------------------------------------------------
// class HybridChoice
//--------------------------------------------------------------------------

VarImplPtr HybridChoise::get() throw (SplitChooser::NoSplitableVar)
{
	VariableSet *sv = m_chooser->getVariables();
	VarImplPtr time = *(m_chooser->getVariables(2)->begin());
	VarImplSet::const_iterator imv, itend;

	if (time->isSplittable()) {
		//m_var = time;
		return time;
	}

	if (m_var.get()) {
		// choice between the variable after m_var and the last one in the model
		imv = m_it;

		for (++imv; imv != sv->end(); ++imv) {
			if (**imv != *time && !(m_chooser->isSolution(**imv))) {
				m_it = imv;
				m_var = *imv;
				return m_var;
			}
		}

		// choice between the first variable in the model and the one before m_var

		itend = m_it;
		++itend;
		for (imv = sv->begin(); imv != itend; ++imv) {
			if (**imv != *time && !(m_chooser->isSolution(**imv))) {
				m_it = imv;
				m_var = *imv;
				return m_var;
			}
		}
	} else {
		// first choice to be made (m_var==NULL)
		for (imv = sv->begin(); imv != sv->end(); ++imv) {
			if (**imv != *time && !(m_chooser->isSolution(**imv))) {
				m_it = imv;
				m_var = *imv;
				return m_var;
			}
		}
	}
	throw SplitChooser::NoSplitableVar();
}

/*
VarImplPtr HybridChoise::get() throw (SplitChooser::NoSplitableVar)
{
	VariableSet *vs = m_chooser->getVariables();
	VarImplPtr time = *(m_chooser->getVariables(2)->begin());
	VarImplSet::const_iterator vis, itEnd;

	if (!m_var.get()) {

		// 最初は時間を選ぶ。
		m_it = vs->begin();
		m_var = time;
		return m_var;
	
	} else {

		// choice from the all variables in the model except time
		vis = m_it;

		for (; vis != vs->end(); ++vis) {
			if (**vis != *time && !(m_chooser->isSolution(**vis))) {
				m_it = vis;
				m_var = *vis;
				return m_var;
			}
		}

		// choice between the first variable in the model and the one before m_var

		itEnd = m_it;
		++itEnd;
		for (vis = vs->begin(); vis != itEnd; ++vis) {
			if (**vis != *time && !(m_chooser->isSolution(**vis))) {
				m_it = vis;
				m_var = *vis;
				return m_var;
			}
		}
	}
	//throw SplitChooser::NoSplitableVar();
	m_var.reset();
	return get();
}
*/
